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We study theoretically vesicle locomotion due to haptotaxis. Haptotaxis is referred to motion 
induced by an adhesion gradient on a substrate. The problem is solved within a local approximation 
where a Rayleigh-type dissipation is adopted. The dynamical model is akin to the Rousse model for 
polymers. A powerful gauge-field invariant formulation is used to solve a dynamical model which 
includes a kind of dissipation due to bond breaking/restoring with the substrate. For a stationary 
i-Q ■ situation where the vesicle acquires a constant drift velocity, we formulate the propulsion problem 

^ | in terms of a nonlinear eigenvalue (the a priori unknown drift velocity) one of Barenblat-Zeldovitch 

type. A counting argument shows that the velocity belongs to a discrete set. For a relatively tense 
vesicle, we provide an analytical expression for the drift velocity as a function of relevant parameters. 
We find good agreement with the full numerical solution. Despite the oversimplification of the model 
it allows the identification of a relevant quantity, namely the adhesion length, which turns out to 
be crucial also in the nonlocal model in the presence of hydrodynamics, a situation on which we 
have recently reported [I. Cantat, and C. Misbah, Phys. Rev. Lett. 83, 235 (1999)] and which 
^2 ■ constitutes the subject of a forthcoming extensive study. 
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PACS numbers 87.22.Bt, 87.45 -k, 47.55 -Dz 



^fj . I. INTRODUCTION 

a , 

Phospholipids vesicles constitute a simple model of cytoplasmic membranes of real cells. A simple model due to 
Hclfrich |Q based on curvature energy has accounted for a variety of equilibrium shapes. The model is based on a 
minimal energy principle Q . Some of the shapes (the so-called discocytes) bear strong resemblance with that of an 
erythrocyte (the red blood cell). Additionally, analysis of flickering (temporal small fluctuations around a given shape) 
of an erythrocyte by Brochard and Lennon @] has been quite successfully described by the Helfrich model including 
hydrodynamics dissipation. The vesicle model has seemed then as a natural candidate, at least in a first stage, for 
dealing with more complex entities such as those encountered in the realm of biology. In that context, however, most 
£NJ i of the features are of nonequilibrium dissipative nature. Very recently several theoretical P-[lO| and experimental pd| ] 
' investigations have been directed along that line. 

We are interested here in vesicle migration, a question on which we have given recently a brief account ||f7|. Despite 
the very complex biochemical behaviour of a cell, cells may also exhibit behaviours where simple physica concepts 
may be evoked. It is well documented that, for example, the migration of the pronephric duct cells in salamanders 
is regulated by haptotaxis. Haptotaxis is a terminology that is used to express the following fact: when adhesive 
molecules are present in increasing amounts along an extracellular matrix (or simply on a substrate in vitro), a 
cell that was constantly making and breaking adhesion with such a molecule would migrate from a region of low 
concentration to an area where that adhesive molecule was more highly concentrated |l^,[l3| . There are also evidences 
that cell migration during embryo development may be guided by an adhesion gradient. In other words cell migration 
is here guided by a purely external physical factor, while the internal structure (the cytoskeletton) is quite unaffected 
, on the time scale of interest. This feature drastically differs from that of a cell belonging to the immune system 
where the cytoskeletton plays a decisive role 14 . Despite the fact that the cytoskeletton in pure haptotaxis does 
not undergo a structural change as is the case during cell crawling, the problem remains very much involved since 
the cell cytoskeletton dissipation may come to the fore as well as an intricate bond breaking and restoring with the 
substrate. We shall consider here a pure vesicle moving in haptotaxis. Our belief is that advancement in this field 
can be achieved only by the progressive refinement of concepts. 

We consider a vesicle moving along the substrate thanks to an adhesion gradient. As the vesicle moves, it generates 
hydrodynamics flow both inside and outside. Hydrodynamics induces nonlocal interactions leading to an effective 
coupling of two distinct regions on the vesicle. In addition, the two monolayers that form the phospholipidic membrane 
might slide one relative to the other. Finally during motion the vesicle forms new bonds ahead and destroys others 
behind, and this process of bond breaking and restoring may be so slow that it may dominate dynamics (see later) . 

This paper should be regarded as using a very simplistic view in the hope of introducing the concepts of migration 
and to exhibit in a transparent fashion the way the problem is addressed. We shall keep the description as simple 
as possible. That is to say: (i) we ignore nonlocality due to hydrodynamics -incorporation of hydrodynamics was 
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briefly discussed in and will be the subject of a forthcoming paper-, (ii) we confine ourselves to a 2D geometry. 
Some kind of dissipation due to bond breaking and restoring is introduced in our model. The adoption of a local 
model (no hydrodynamics) allows one to quite easily obtain analytical results and thus to extract some key ingredients 
about migration -especially the role of the adhesion area (length in 2D)- which turns out to be crucial also when 
hydrodynamics is included. 

The scheme of this paper is organized as follows. In Section II we write down the equations of motion and comment 
them. In section III we present a forward time integration and present the main result. Section IV presents the solution 
of the stationary system in a form of a nonlinear eigenvalue problem, where the drift velocity is the eigenvalue. In 
section V we give an analytical solution. A conclusion and a discussion is presented in section VI. 



II. EQUATION OF MOTION 



A. Parameterization 



We consider an adhering vesicle, deposited on a flat substratum which is oriented by its normal vector y (Fig.l). 
The x-axis is along the wall and represents the direction of vesicle motion occurring by convention from left to right. 
As stated before we confine ourselves to two dimensions. That is to say, the vesicle morphology is invariant in the 
z-direction, similar to a tubular vesicle. 

The interaction between the vesicle and the substrate is taken into account by introducing an adhesion potential. 
The range of the potential in realistic situations (typically several nm) is small in comparison to the vesicle size 
(several fim), so that it is justified in practice to consider a contact potential, unless otherwise indicated (see later). 
The energy interaction is then zero if y > and is non vanishing only close to contact (if y = 0). At the junction 
point between the free part of the vesicle (whose length is denoted as L*) and the adhered part (with length L a dh), 
the potential undergoes an abrupt change. The contact between the membrane and substrate occurs at two well- 
defined points xi at the left and x-i at the right. These parameters are related to the total length of the curve L by 
L = L a dh + L* = (x2 — xi) + L* . We use an intrinsic representation of the vesicle contour by introducing tp(s), the 
angle between the outward normal and the y-axis, and s the arc length, as shown on Fig. 1. We only need to consider 
the function tp(s, t) from s — to s = L* corresponding to the contact points X\ and X2, respectively. Because of the 
contact potential character the adhesion length is completely fixed if the two contact points x\ and xi are known. 
Thus the vesicle shape and its dynamical properties (like the propulsion velocity) are known if x\ , X2 and the function 
ij}(s,t) are determined. The demand that the parametrization of the vesicle be compatible with the adhesion on the 
substrate it fulfilled by the two geometrical constraints : (i) the distance between both contact points, Xa — Xi, must 
coincide with the adhesion length L a dh = x% — Xi, (ii) their vertical coordinates y\ and y 2 must have the same value, 
J/2 — J/i = 0. These two constraints can be expressed in terms of ip(s,t). For that purpose we use the relations 

dx , dy 

-=cos<A, ^ (1) 
which allow us to write the two constraints in the following form : 



f dx f 

/ — ds = / cos i/)(s) ds = x 2 - xi = L a dh , (2) 
Jo ds S Jo 



- sinip(s) ds — y 2 - yi — . (3) 



These are the geometrical constraints. In order to describe vesicle dynamics, we need a dynamical equation for the 
evolution of ip(s,t). A movement of the vesicle (due for example to an adhesion gradient) is limited by dissipation 
(such as hydrodynamics etc.). The vesicle reacts to any deviation from equilibrium by its internal forces (bending, 
possible stretching -or resistance to stretching-). Let us first discuss these forces. 



B. Energy and forces 

All the relevant membrane properties are summarized in the following energy, expressed in 2D, with the dimension 
of an energy per unit length : 
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£ = / ft — ^-^ds- / w(ir)dizi + / ((s)ds + p S . (4) 
Jc 2 in Jc 

The first term is the well known Helfrich curvature energy, with the rigidity k, the curvature c = dip/ds and the 
spontaneous curvature c s jjj. Because of the 2D-specific conservation law J cds = 2ir, any curve displacement leaves 
unchanged the energy terms associated to the spontaneous curvature. We can thus omit the term associated with 
c s . The second term expresses the adhesion energy and is only integrated on the adhered part of the curve. As we 
are concerned with an inhomogeneous substratum, the contact potential depends on the variable x and is denoted by 
— w(x) (with w > 0, meaning that adhesion is favorable). Finally the last two terms ensure the length and surface 
conservation, respectively. The membrane is a two dimensional incompressible fluid. The phospholipid exchanges 
with the solvent is virtually absent, and the area per molecule on the vesicle remains constant. This leads to the local 
length conservation (in the 2D language). The variable C is a local Lagrange multiplier which enforces the arc length 
ds to a constant value. The enclosed surface S conservation is a consequence of the membrane impermeability and 
fluid incompressibility. It is ensured by the global Lagrange multiplier p. The interpretation of £ and p as a tension 
and a pressure will be discussed latter. 

The functional derivative of the energy (Q) provides us with the various forces acting on the membrane. 

(i) Curvature forces 

Under the assumption that the membrane is completely flat on the adhered part, we reduce the integration domain of 
the first energetic term in eq.(^) to [0, L*]. Making use of the relations t = r' and n = — (l/c)r" (the prime designates 
derivative with respect to s) which are the membrane tangential and normal unit vectors, we obtain for the curvature 
energy E c 

k f L ' fd 2 r^ 2 



When taking the functional derivative of E c with respect to the position, care must be taken. Indeed the arc length 
element must also undergo a variation. A convenient formulation avoiding confusion rests on the introduction of a 
general parametrization a £ [0,1], related to s by the metric g = (ds/da) 2 , and is time-independent. The energy 
expression becomes fqjlq] 
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d 2 r\ 2 1 / d 2 r dr x 



da 2 J q \ da 2 da 



- z ' 2 da . (6) 



The functional derivative, though straightforward, may be too lengthy if one does not take care in regrouping ade- 
quately various terms as explained in |8|. The result can be written in a simple form: 



5E C fd 2 c 



Sr(s) 



(7) 



The curvature force is, as expected, free of any tangential contribution. The first term in eq. ([?]), involving the 
second derivative of the curvature, tends to keep curvature repartition as homogeneous as possible. It is also present 
in 3D under the more complicated form of the Laplace-Beltrami operator. The second term proportional to c 3 is in 
contrast 2D-specific. It tends to increase the size of any convex shape. Note that in 3D the curvature energy is scale 
invariant, which implies a vanishing curvature force of this type on a sphere. The difference between 2D and 3D can 
be explained in the following way. Let us consider a finite cylinder of length H and radius R <C H, which constitutes 
a good approximation for an infinitely long cylinder. In order to make the cylinder "closer" to a sphere, which is the 
corresponding 3D equilibrium shape, the curvature force would tend to increase the radius and decrease the length so 
as to bring the cylinder shape as close as possible to a sphere. This gives an intuitive picture of the c 3 term in 2D. 

In the discussion above we did omit the boundary contribution when taking the functional derivative. Since this 
point is a bit subtle we have postponed it to the end of this section. 

(ii) Length and surface constraints 

On the free part of the curve, the force which is associated to the first Lagrange multiplier £ is obtained upon functional 
derivation of E\ = f (ds. The result is 

{ l = -6E l /6r(s) = -c((s)n+^t. (8) 

ds 
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The normal component is easily identified as a Laplace pressure, whereas the tangential one looks like a Marangoni 
force (which is encountered when surface tension is inhomogeneous). However, £ is not exactly similar to a surface 
tension as for an interface between two fluids. The "tension" £ is not an intrinsic property of the membrane. It adapts 
itself to the other forces in order to maintain the local length fixed. In other words, the problem is implicitly written 
in a thcrmodynamical ensemble with fixed length. This differs from the usual problem for fluid or solid surfaces where 
the surface tension is fixed instead. Thus £ is a variable that must be determined self consistently as a Lagrange 
multiplier, by use of the constraint equation (see appendix in Ref. JTq|) : 



= 



d(ds) f dv t 



= I -7T- + CV n I ds 

OS 



(9) 



This relation (g|) simply expresses the condition of vanishing velocity divergence on the curved contour of the vesicle, 
which is precisely the incompressibility condition in the 2D fluid constituting the membrane (written here in one 
dimension). A more intuitive way of viewing expression (0) is presented on Fig- 2. The Marangoni term is the only 
tangential term among all membrane forces (see eqs. (^,jsljlO|) . It is seen from (||) that the Lagrange multiplier must 
be uniform at equilibrium. For sake of simplicity and in order to get more insight into analytical understanding, a 
uniform value will be assigned to Cj even out-of-equilibrium. A discussion of this point will be presented in section 
VI. This assumption imp lies some consequences on dynamics (and especially on the tangential velocity) which will be 



presented in section [I C 



Finally we have to consider the force associated to E s = pS 

f s = -SE S /Sr{s) = 



-p n 



(10) 



The Lagrange multiplier p depends only on time; it enforces a constant area. Two physical interpretations can be 
invoked depending on the situation under consideration. Either we consider an impermeable membrane, and p would 
be the hydrostatic pressure difference between outside and inside; or we choose a model of permeable membrane and 
p would play the role of an osmotic pressure. Both models are equivalent as long as we do not consider hydrodynamic 
flows. 

(iii) Adhesion forces and boundary terms 

The functional derivative induces boundary terms at each contact point. The additional variation SE h c and SE h w l for 
the curvature, adhesion and tension energies, associated with a small displacement Sr of the contact points is given 
by (see §jl|) 



8El 



Sr- 



Sr ■ 



V9 

Sr ■ Ut + {(- w(x))x 



dc 

OS 



L* 



(11) 

(12) 



Following the definition of these boundary points, they remain on the substrate. Thus, the accessible values for Sr 
is then reduced to Sr oc x. Additionally, in order to keep the curvature energy finite, we impose a vanishing value 
for the contact angle <j> between the membrane and the substrate (see Fig. 3). Within our formulation, this constraint 
does not follow from the energy minimization and has thus to be added into the physical model. More precisely, at 
the discontinuity point (say X2) one has to add to the Helfrich energy a term of the form K(Aip / As) 2 which informs 
us on how would the vesicle on the adhered part feels, so to speak, the behavior of the vesicle at the junction point 
on the right side. Across the contact point of a vanishing extent, As — >■ 0, while the angle, if it had to have another 
value than zero, would make a jump leading to an abnormally increasing curvature energy. We must then impose 
a vanishing contact angle. These various conditions (motion along the wall and a vanishing contact angle) lead to 
n = — y, t = —x and Sr oc x. It follows then that the term proportional to Sr in eq.(O) vanishes automatically. The 



second term becomes k/2 \c\Sx2 — c\8x\j with c\ and C2 the curvatures at the left and right contact points. These 
terms are counterbalanced by adhesion and tension terms (eq.(O)) leading to the relation 



(13) 



T U__^) |X 



where the — and + signs refer to the rear and fore contact points represented by the subscript z = 1,2. At equilibrium, 
we recover here the relation c = y /, 2w/n jlT]] . 

The energy variation given by eq.(13) can not really be identified as a physical force. It corresponds indeed to a 
geometrical point displacement. The "force" orientation is here parallel to the substrate, whereas the real force acting 
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on the contact point, considered as a material point, is expected to be normal to the substrate. As we have seen above 
the curvature forces are indeed normal when applied to a an adjacent piece of the membrane (see eq. 0). The present 
" force" has the meaning of how much energy would be involved in displacing the contact point from one position to 
another. That geometrical point is by its very nature sitting on the substrate, so that the "force" associated with its 
displacement is naturally tangential. 

We find it worthwhile to make a short digression. Suppose that the angle is not fixed to zero as we did above. More 
precisely suppose that the rigidity is so small or the adhesion is so large (see below what does this mean) then the 
vesicle will be so tense that it would look like a droplet outside some length scale I to be determined below (of course 
within that scale, which is sufficiently close to the substrate, the matching must be tangential). If we do not assume 
a value for the angle (that is no relation between n or t with x and y), and set <5r cx x we find from (|l3|) that c = at 
the contact (which means a straight line at the contact) and that the angle between that line and the substrate obeys 

Wi = (1 - cos((j> % )) C 

<f>i ~ \j — ^— ( for small angles ) (14) 

which is nothing but the Young condition. We have neglected ndc/ds in comparison to w. The justification is as 
follows, ndc/ds ~ kcq/£, where cq is the true contact curvature given by y/2w / n. The approximation is legitimate 
provided that the length scale i ^> yj k/w. The length \J n/w is the radius of curvature at contact. If the scale of 
interest is outside that internal region, then the droplet limit is justified. It must be emphasized however that the 
effective contact angle is not an intrinsic property of the adhered membrane, as for a droplet, but it is linked to other 
parameters (rigidity, the vesicle scale-on which depends £— , etc... ). In particular, the tension £ is fixed by the reduced 
volume, which is a global property of the vesicle : different vesicles of the same phospholipid composition, but with 
different sizes, may have different contact angle on the same substrate. 



C. Dynamical equation 

An important point which must be emphazised when dealing with dynamics is the identification of the dissipation 
sources. These are the following: (i) the dissipation in the membrane via molecule rotations (very much like liquid 
crystals where dissipation is characterized by the Leslie coefficient), (ii) hydrodynamics flows inside and outside 
the vesicle, (iii) friction between the monolayers, and (iv) bond breaking and restoring with the substrate. It is 
well known that dissipation associated with rotation (internal dissipation) is negligible in practice and for free 
vesicles (no substrate) hydrodynamics seems to be the most important dissipation. Hydrodynamics induces nonlocal 
interactions J?J and this will be dealt with extensively in a forthcoming paper. Our wish in this paper is to present 
a pedestrian model, namely a local one, which allows for a complete analytical solution that will help to identify 
some key ingredients in the migration process. A specific dissipation with the substrate will be introduced later. 
For the moment we confine our description to the free vesicle case. The local model to be presented here is similar 
to the so-called Rousse model in the community of polymers. Indeed, for a one dimensional contour in a three 
dimensional space dynamics becomes local even in the presence of hydrodynamics |l9[ . 

The best way to introduce the dynamical law is to consider a dissipation function, proportional at each point to 
the square of the velocity : 



F <l = jJ . (15) 

The coefficient r/ is here an effective viscosity and has the dimension of a viscosity per unit length. Its numerical value 
is estimated by rj = r) wa t/R ~ I0 2 kgm~ 2 s~ 1 , with rj wat the water viscosity and R a typical vesicle size. 
Neglecting inertial terms, the Euler-Lagrange equations become then 

- — = — => 7?V = f . (16) 

or ov 

As expected, we find a local proportionality between the membrane velocity v and the membrane force f , which is a 
nonlinear function of position. In the present picture where the effective tension £ is space-independent no tangential 
force appears so that physics will only fix the normal velocity, while the tangential velocity has no physical meaning 
as described below. 

(i) Normal velocity 

The normal membrane force is (see eq. (0, |[ |o|)) : 
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d 2 c c 3 



From the dynamical law ( |16| ) and the membrane forces expression (Jl7|) we obtain the normal velocity as a function of 

i>{8) : 



Vn(s) = fn = ~ 



d 3 t/j l(d±\ _(d±_P 
ds 3 2 \ ds J k ds k 



(18) 



It is convenient to write the dynamical equation in terms of the angle tp and not the curvature c. The reason is that 
the boundary conditions are written naturally as a function of tp (tangential matching, tp = ±7r, and contact curvature 
dtp/ds = y/5w/K). 

(ii) Tangential velocity 

There is only one tangential contribution to the membrane forces, d(^/ds, which is zero with the assumption of a 
uniform tension (see eq. |J). This implies that only the total length is conserved, and not the local one. A dilatation 
of a part of the membrane is then permitted, as long as the remaining part of the vesicle is contracted in order to 
keep the total length unchanged (see Fig. 4). Within this approximation, there is no energy variation associated to 
tangential motion, and therefore no forces. In other words we consider the vesicle contour as a mathematical curve, 
loosing the concept of density : only the shape matters, independently of the points distribution on the curve. 

We could equivalently assume that only the normal velocity contributes to dissipation. In that case the dissipation 
function ( |l5| ) would take the form Fd n — "n/2 J \v n \ 2 ds and the equation of motion ( |l6| ) becomes f = r)v n ii. The 
tangential force d(^/ds must then vanish and we get automatically that £ ^constant. Thus our assumption of a global 
Lagrange multiplier can also be viewed as the result of a dissipation due uniquely to normal displacements. 

As we have already mentioned, tangential displacements do not induce a geometrical change. If the tangential 
velocity has no physical meaning (as is the case with a constant Q its choice should not affect the physics. We are thus 
at liberty to choose one which is convenient (very much like a gauge- field invariant formulation in electrodynamics). 
The choice of a gauge is interpreted as a reparametrization of the curve. As seen below the tangential velocity is fixed 
by the normal velocity once the gauge is specified, but there is naturally no feed back of that tangential velocity on 
the normal one (the physical one). This is the crucial difference between this non physical tangential velocity and 
a physical one that may arise in the general case (as discussed in a forthcoming paper). It must be noted however, 
that whether a tangential velocity is of physical nature or not, the knowledge of the normal velocity is sufficient to 
describe vesicle dynamics. It is thus only via its influence on the normal velocity that a physical tangential velocity 
would affect the physics (see below). Still in that case we can introduce a second tangential velocity of geometrical 
nature that corresponds to the displacement of the representative points of the curve and not to the material ones 
which are affected by the physical tangential velocity. 

In the present model the most convenient parametrization requires a homogeneous points distribution along the 
free part of the curve, which is expressed as d/dt (s(a)/L*) = 0. This provides the expression for the "non physical 
tangential velocity" (see appendix in Ref. p6|) : 

v t (s) =«t(0) -jT cv n ds + ^ ^jf cv n ds + v t (L*)-v t (0)J . (19) 

If the free length L* remains constant during the motion, as happens for a stationary regime, eq.([l^) fixing the gauge 
imposes nothing but a constant distance between two consecutive points on the vesicle. The local length conservation 
(^), which is physical, seems then to be implied by a gauge!. In reality, once we have adopted a contant tension 
-implying a vanishing tangential physical velocity-, any point distribution is of purely geometrical nature, and we 
could impose another gauge than the above one, without affecting the physics; the above tangential velocity does 
not act on the normal velocity. Had we considereed £ to be non contant, we would then have obtained a physical 
tangential velocity, which would act on the normal one; use of (^J) fixes £(s) which in turn acts on the value of v n , and 
then on physics. In the simplistic model we adopt, the tangential velocity is determined a posteriori, independently 
of the normal velocity. That is why it is only a non physical reparametrization, a "gauge". A remark is in order : 
in this situation the question of a rolling or sliding motion does not make sense, since both motions differ only by a 
tangential velocity. 

The membrane velocity is given as a function of tp(a), a being the auxiliary parametrization of the free part of the 
vesicle, running from one contact point to the other. In order to obtain a closed system we need a relation between 
the evolution of tp(a) and the velocities. The temporal derivative of ip, for a given a, is presented in the appendix of 
Ref. M 
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-dt) a = CVt -^i- (20) 

The last step is the determination of the boundary conditions at the contact with the substrate. 

(iii) Contact points velocity 

The motion of the contact point is governed by a binding/unbinding mechanism, implying a dissipation law that 
differs from the bulk dissipation. The most natural way for introducing a dissipation law is the following (with T a 
phcnomenological dissipation coefficient) 

Using the energy variation ([L3j) we get the following dynamical law, with Wi — w(xi) and Vi — dxi/dt 



C1 = V K ' C2 = V K ' (22) 

which constitute the dynamical boundary conditions. These out of equilibrium values for the curvature are quite 
intuitive : the unbinding delay at the rear point induces a larger curvature than at equilibrium, whereas the binding 
delay at the fore point induces a smaller curvature. 



III. TRANSIENT BEHAVIOR 



The formalism presented in the first part lends itself very well to analytical computation and stationary shape 
determination, as will appear in the following paragraphs. Nevertheless, having access to the transient process is 
highly desirable. In particularly, it checks the dynamical stability of an eventual stationary behavior, obtained after 
a relaxation. The successive vesicle profiles are determined by a direct numerical implementation of the dynamical 
equations ([l8|), (|l9|) and (^0|). Unfortunately, numerical instabilities are difficult to avoid around each contact point 
(due to a contact adhesion potential). A smoother model, without discontinuities, is more convenient for such an 
approach. For this reason, in this paragraph devoted to transient processes, the adhesion potential will be supposed 
to be of small, but non vanishing, range. We rapidly summarize below the small technical changes arising from this 
model modification. The chosen potential profile is 

W (r)=w (l + u x)(4-^-), (23) 

y y 

with the new length yo fixing the characteristic distance between the substrate and the membrane, and w(x) — 
— Wo(l + uqx) the minimum of the potential interaction, occurring for y = y$. It plays the role of the previous w(x). 
It depends linearly on x with an adhesion gradient uq. The distance yo is chosen of the order of lOnm, which is 
small enough in regard of the vesicle size to introduce only small variations between both models. In this case the 
parametrization is performed on the whole closed curve, and the boundary terms (eqs. |22] ) are no more relevant. The 
adhesion forces f w are obtained by functional derivation of / w(r)ds leading to 

f w = -(cw + Vw ■ n)n . (24) 

Additionally the gauge condition fixing the tangential velocity eq. ( jl9| ) is simplified : the velocity v t (0) is supposed 
to be zero, without loss of generality, so the first term disappears ; the last term is proportional to the length variation 
of the total parametrized curve, which is zero because we consider the complete profile and no more the free part of 
the curve. Thus we obtain for the gauge, replacing eq. (|ll]) : 

v t = - / dsv n c . (25) 



Using equations (24) and (p5|) we finally get the dynamical equation for r 



dr(a,t) 

Ft 



d 2 c 1 

k ( — 7? H — c 6 I — cw — (Vw • n) — p — (c 
1 ds 1 2 ' 



n + v t t . (26) 



The Lagrangian multipliers are determined from the following constraint equations : 
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dL f 

— = / cv n ds = , (27) 
^ = / M« - . (28) 



pressure and tension va. 



The normal component of the velocity in eq. (26) will be denoted by convention as v n — w„ — p — cC- With this 
notation the equation (2^) appears as a very simple linear equation system in £ and p. Its solution provides the 
ues : 



C= Ml=i£>M> , (29) 



P = -CM + (O . (30) 

with the average defined by 



f 



L 

ds---. (31) 



We are now in a position to deal with the numerical anlysis. The dynamics is overdamped and for this reason 
local in time. Starting from an arbitrary profile, forward time integration provides us with the vesicle evolution. We 
have checked that (i) for a free vesicle (no substrate) the shape (with no external force) tends towards that obtained 
by direct energy minimization, (ii) we have also checked that for a homogeneous substrate an arbitrary initial shape 
evolves after some time to the shapes obtained in [ p"7| by direct energy minimization. 

Let us now turn to the non-equilibrium situation ensured by an adhesion gradient. Starting from an initial shape, 
the vesicle acquires a non-symmetric shape and moves in the gradient direction. After transients have decayed the 
vesicle acquires a permanent regime with a constant velocity. Figure 5 shows the shape evolution. 

IV. STATIONARY MOTION : DIRECT NUMERICAL SOLUTION 

The formulation of our problem in terms of a direct stationary problem is very convenient both for a systematic 
study of the velocity evolution as a function of various parameters. It will also allow us to present a simple analytical 
solution. It is convenient here to come back to the contact potential model. For a vesicle which has attained a 
stationary shape and velocity V the equations become steady with V as an unknown parameter. 

For a stationary motion along the x-axis, normal and tangential velocities can be written as functions of the angle 
if) and of the translational velocity V : 

v n = Vx ■ n = V sin ip , (32) 
v t = Vx ■ t = V cos ip . (33) 

The shape and velocity are entirely determined from the relation between normal velocities and forces. The equation 
of motion is obtained from eqs. dl8,B2|): 



V sin ip = — 



d 3 tp 1 fdip\ 3 (dip p 
ds 3 2 \ ds J k ds k 



(34) 



Let us present briefly a counting argument showing that the problem is well defined. Equation ( p4|) is a nonlinear 
third order differential equation for if), with 3 parameters to be determined : p and V. So we need 6 "informations". 

We have the following equations at our disposal : 
• 2 geometrical constraints (eqs. (||) and (||)) : 

cos if>(s) ds = L a dh , / — sinip(s) ds = (35) 



• 4 boundary equations corresponding to the contact angles and their first derivatives (the dynamical contact curva- 
tures ci and C2 given by equation (p2|)) 

dip \ dip \ 

1p! = 1p{S = 0) = -7T, 1p2 = 1p(s = L*) = 7T , —J = Ci , J = C 2 (36) 
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• 1 equation ensuring that the enclosed surface is equal to the prescribed area. 

• 1 equation ensuring that the total length of the curve is precisely the prescribed one, L, which is related to the two 
other lengths by 

L = L* + L adh . (37) 



There are thus 8 conditions, for only 6 informations needed. The system seems then to be overdctcrmincd. This is 
not the case. Indeed it must be noted that the problem involves additional unknowns which are L* and L a dh- So in 
reality we have 8 unknown parameters as well. The problem is thus well defined. 

Once the shape is determined we must in principle evaluate the area and change the parameter p until the area 
coincides with the prescribed one. But since the area is a conjugate variable to p we can fix p -which is more 
convenient- and this will fix some area that is treated as free (not imposed in advance) . Additionally we are at liberty 
to prescribe L (that fixes some length scale). L a dh can then be determined if L* is known; L ac ih can thus be removed 
from the problem upon using eq.(|37|). The first constraint ( |35| ) becomes then 

cosip(s) ds = L — L* . (38) 

In other words prescribing the total length to L and the pressure to p lowers the number of unknowns by two. This 
is so because we do not want to have a specific area, and that L* and L ac ih are not independent if we treat the total 
length as known. That is to say we have finally 6 fixed boundary conditions or constraints ( |35f|36| ) and six parameters 
which are L* , V and £, plus three constants of integration due to the third order differential equation (|34). 

A convenient way to solve a differential equation of order n is to transform it into a set of n first order coupled 
differential equations. For that purpose we set fx = ip, fx = fa and /2 = fa where the dot stands for d/ds. Equation 
( |34| ) then provides us with the expression for fa = d 3 ?p/ds 3 as a function of fx, fi : 

h = »?/« P2 sin h ~ / 2 3 /2 - p/k + hPzl* = F . (39) 

In order to make visible the quantities which are treated as unknown parameters we shall use the symbols Pi (with 
i = 1,2...). As stated above there are three parameters L* = P%, V — Pi and ( — P3. Solution of a set of three 
equations involves three integration factors. This means that we have 6 unknowns, as argued in the last paragraph. 



Four physical conditions at the two end points (see eq. ( pq)) are known. Two constraints are imposed (eq. (35)), and 
this makes the problem well posed. Note that conditions (|35|) have an integral form. We find it convenient to rewrite 
them in a differential form. It is easy to realize that by setting 

fi(s)= [ siniP(s')ds' , / 5 (s) = f cosip(s')ds' , (40) 
Jo Jo 

we can write 

fx = sin?/; = sin/i , f 5 = cost/; = cos/x . (41) 

These two functions obviously obey /4(C)) = 0, /5(C)) = , whereas at the second boundary we must impose fx(L*) = 
and fz{L*) = L — L* in order to fulfill the two constraints ( ^5| , [j8| ). This trick is performed at the expense of two 
additional functions fi and f§ (whose determinations involve two integration constants). We have thus augmented 
our system by 2 differential equations of first order. The two additional integration constants are precisely fixed by 
the demand /4(D) = 0, /s(0) = 0, whereas the conditions fi(L*) = and f${L*) = L - L* are substituted to (|3^). 
Finally the shooting NAG code used here requires to invoke the boundary conditions for each function /j, with i < 5. 
The boundary conditions for each function is invoked above, except for /a which represents the second derivative of 
ip. This quantity is not known at the boundaries and there is no constraint to be imposed on it. Let P4 and P5 
denote the values of f^ at the two end points. We can thus invoke the boundary conditions of f^ and the boundary 
values are quantities which are to be determined. That is to say we introduce two conditions with two additional 
unknown parameters. We have then in total ten unknowns and ten conditions. Cast into this form our formulation 
can straightforwardly be implemented into a NAG code (code D02HBF). 

In summary the problem to be solved can be written in a standard boundary value problem with unknown param- 
eters : 



(42) 



fx 


= h , 


/i(0) 




h(pl) 


= TT 


h 


= h 7 


/ 2 (0) 


= ci , 


/20P1) 


= c 2 


h 


= F, 


/s(0) 




fs(Pi) 


= P 5 


h 


= sin fx , 


/ 4 (0) 




U(Pi) 


= 


h 


= cos fx , 


/b(0) 


= 0, 


h(pl) 


= L-Px 
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Once the problem is solved the vesicle shape is obtained by making use of equations ([I]) . 

If wi = u>2 j the vesicle is at equilibrium on a homogeneous substrate and one obviously expects a vanishing velocity. 
This comes out automatically from the above formulation. If we were interested from the beginning in an equilibrium 
problem, we would then not have introduced V as an unknown parameter. In that case because the profile is symmetric 
the second condition (|3^) is automatically satisfied, and we would then be left with with nine conditions for nine 
unknowns. 

When wi u>2, there is no equilibrium solution for the vesicle, which has to move towards the stronger adhesion 
region. If we impose a vanishing velocity there is no way to fulfill the second condition ( |35|) (a typical profile would 
be the one shown on Fig. 7) where starting from one end we arrive at the other end at a different height. Arriving at 
the same height can be achieved only for a specific velocity (or at most a discrete set of solutions), the one we are 
seeking. Thus the second condition of eq.(|35|) (which is parametrized by the set of Pi) can be viewed as 'quantization' 
condition. This is a nonlinear eigenvalue problem of Barenblat-Zcldovitch type. 

The numerical solution reveals an out-of-equilibrium shape which is significantly different from the equilibrium one, 
as shown on Fig. 6. We note that the curvature in front of the vesicle is higher than the one behind. The reason is 
that the adhesion energy is higher in the front part, so that the curvature/adhesion balance allows a higher curvature 
(the vesicle looses curvature energy at the expense of a stronger adhesion) . 



V. STATIONARY MOTION : AN ANALYTICAL SOLUTION 



The advantage offered by the simplistic picture of our model is the possibility to provide analytical results and thus 
to shed light on the physical processes that are involved in the problem of vesicle propulsion. It turns out that the 
equation of motion (|34|), if multiplied by d 2 ip/ds 2 , possesses practically a first integral 



V 



L' 



d 2 1p 



ds 2 



V 



C, 



^ sm(iP)ds = ~ ( [(dc/d S ) 2 /2]\ + ±[c 4 /4] 2 - ~[c 2 /2]l - V -\c 



(43) 



Each r.h.s. term has an explicit form as a function of the contact curvatures (which are known), and of their first 
derivatives, which have only negligible contribution for swelled vesicle. The l.h.s. term can be evaluated for a vesicle 
shape close to a circle. The calculation is sent into appendix and leads to 



d 2 ip 
~ds 



2 sin(ip)ds 



4n 2 



. L, 



dli 



(44) 



Using the dynamical values for the contact curvature given by eq. (|22|) , we obtain an explicit expression for the velocity 

~2l 



V = 



(2 



nrvLadh V 2 K 



(45) 



In the simple case where T — (no dissipation associated with the substrate), expression (|4^) becomes explicit and 
provides a good agreement with numerical solution (see Fig. 8). The analytical expression for the velocity involves 
only known parameters, except L a dh and C F° r the comparison between numerical and analytical results, we took 
their numerical values. 

Limit <5w <C 1 

Another interesting limiting case concerns the small adhesion difference. Expansion of the numerator in cq.(|45|) to 
leading order in Sw yields 



V 



Sw 



rj/A + T 



A = 



w R 2 



dh 



1 - 



P 



(46) 



where R — L/2n. The influence of the two dissipation coefficients appears then clearly. It depends on the quantity 
A, proportional to the ratio R/L a dh- The bulk dissipation increases with the adhesion length, which seems to be 
a very robust result, as encountered in the model including hydrodynamics dissipation (ic| . The local dissipation 
represented by T does not depend on the adhesion length. Only the two contact points matter. Note also that the 
effective dissipation is r]/A + T. The bulk dissipation r\ /A and the contact one T play a role of resistances (in an 
electric analogy) which would be mounted in series. 
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VI. DISCUSSION AND CONCLUSION 



This paper has given a first extensive presentation of the problem of vesicle migration in haptotaxis. We have 
reduced as much as possible the complexity of the problem in order to gain some analytical approximate results. For 
that purpose we have neglected hydrodynamics which induce nonlocal interactions, and adopted a local model of the 
Rousse type. The full dynamical problem has been solved by adopting a powerful gauge-held invariant formulation. 
The dynamical code could account for the transient and the evolution towards a steady-state solution. In that 
context an introduction of an adhesion potential with a finite, albeit small, range has proven to be necessary in order 
to circumvent numerical instabilities related to the motion of the contact point. This code has the advantage of dealing 
with various problems not leading necessarily to permanent motions. For a stationary situation we could cast the 
problem into a standard boundary value one where the migration velocity appeared as an eigenvalue. This problem 
is akin to the nonlinear eigenvalue problem of Barenblat-Zeldovitch type. A counting argument showed us that the 
velocity should belong to a discrete set, only one of them has been identified; we speculate that the solution is unique. 
The problem could be systematically solved in a fully intrinsic representation of the contour. For a rather tense vesicle 
we have provided an analytical solution which is in a good agreement with the numerical one. We have identified the 
role played by the adhesion length in selecting the magnitude of the migration velocity even if no dissipation with the 
substrate is included. We have also shown that the bond breaking/restoring dissipation and the (effective) bulk one 
are additive in a way which is analogous to the problem of electrical resistances in series. Bulk dissipation dominates 
when the ratio of the bulk dissipation coefficient to the contact one exceeds a certain limit, which depends in an 
intricate manner on various parameters. For real situations, vesicles, and cells in general, are suspended in aqueous 
solutions. It is therefore highly important to include hydrodynamics. Moreover the Lagrange multiplier £ is a local 
quantity. We have recently given a brief account on these questions fic| l . An extensive discussion will be presented in 
the near future. 



APPENDIX A: DERIVATION OF EQUATION 44 




We write ip = 2tts/ L* — 7r + e, which implies to first order in e' 




We then make use of the following relation: 




The integral between and L* of the l.h.s. trem vanishes. We then obtain 




The sought after relation has then the form 
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Fig. 1 Notations used in the text. 

Fig. 2 A geometrical explanation of the arc length variation with time. 
Fig. 3 Force equilibrium at the fore contact point in the small rigidity limit. 

Fig. 4 Translation of a circle obtained with a purely normal motion. The right part is dilated, whereas the left part 
is contracted. 

Fig. 5 Successive vesicle profiles. The first one with open circles is an arbitrarily chosen initial shape. It relaxes to a 
permanent shape marked by filled circles on a inhomogeneous substrate. 

Fig. 6 Out of equilibrium adhering vesicle profiles. V is measured in units of 100/zm and W in units of 10~ 4 mJ/m 2 . 
Fig. 7 Geometrical constraint on the curve. 

Fig. 8 Evolution of the vesicle velocity as a function of the adhesion difference. 
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L= Ladh+L 





S (surface) 



Ladh 



low adhesion 



great adhesion C antat.Misbah 




Cantat.Misbah 



AB = ds = a R 

A"B" = a (R + Vn dt) = ds (1+ Vn dt / c) 
A'B'=A"B" + (Vt(B)-Vt(A)) dt 

= ds (1+ Vn dt/c) + ds d(Vt)/ds dt 



d(ds)/dt = ds (d(Vt)/ds +c Vn) 



X1,Y1 



Cantat.Misbah 



X2. 



